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SUMMARY 


When a laminated composite Is subjected to compressive loads, a 
delaminated region may buckle. This causes high interlaminar stresses at the 
delamination front and the delamination may grow. The effect of various 
parameters on instability-related delamination growth was studied analyti- 
cally. The configuration studied consisted of a thick composite laminate with 
a single through-width delamination located near one surface. Both mechanical 
and thermal loads were considered. All conclusions were based on the assump- 
tion that Gj and Gjj govern delamination growth. An approximate super- 
position stress analysis was developed which gives closed form expressions 
for Gj and Gjj. The simplicity of the analysis permitted examination of 
numerous configurations. Both G-j- and G^j were found to be very sensitive 
to delamination length and location through the thickness. The magnitude of 
Gj was also very sensitive to initial imperfection. Critical loads for 
delamination growth were calculated based on three growth criteria. Large 
differences in the predictions highlight the need for a verified mixed-mode 
delamination growth criterion. 
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INTRODUCTION 


In laminated composites under compression loads, a delaminated region may 
buckle. When this happens, high interlaminar stresses can exist at the delam- 
ination front and the delamination may grow. Such delamination growth is 
referred to herein as "instability-related delamination growth" or IRDG. The 
mechanism of IRDG has been observed in a variety of specimens. Under compres- 
sion fatigue, specimens with circular holes develop delaminations near the 
holes. These delaminated regions buckle, which causes IRDG (ref. 1). Delami- 
nations in impacted specimens also have been observed to buckle and grow 
(refs. 2-6). In several studies, specimens with single delaminations have 
been made by imbedding Teflon or Kapton film in the laminate during manu- 
facture (refs. 3, 7-13). When tested in compression, the specimens exhibited 
delamination growth caused by local buckling of the delaminated plies. 

The prediction of instability-related delamination growth is complicated 
by the highly nonlinear character of the deformation. At the heart of the 
IRDG mechanism is a postbuckling process. The postbuckling behavior of a 
column or plate is very sensitive to the in-plane dimensions, the flexural and 
in-plane stiffnesses, and the initial imperfections. Hence, the growth of a 
delamination due to local buckling is expected to be sensitive to these same 
parameters. Other parameters expected to affect IRDG are the type of loading 
(i.e., mechanical or thermal), and, of course, the interlaminar fracture 
toughness. 

The configuration studied (fig. 1) consists of a thick composite laminate 
with a single delamination located near one surface. The delamination is rec- 
tangular and extends across the entire specimen width. This configuration is 
perhaps the simplest that exhibits IRDG. Even so, the response to mechanical 
and thermal loads is still very complicated. Reference 14 presents a closed 
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form analysis of this configuration which calculates total strain-energy 
release rate. But studies (e.g., refs. 10, 15) have shown that the mode I 
fracture toughness is much lower than the mode II fracture toughness. Hence, 
an analysis is needed to calculate the components of strain energy release 
rate (Gj and Gj j) . References 7 and 8 used geometrically nonlinear finite 
element analyses to calculate Gj and Gjj for several configurations. But 
the insight provided by strictly numerical analyses is very limited, since the 
results apply only to the particular configuration analyzed. Reference 8 
dissected the original nonlinear configuration into linear configurations with 
nonlinear ly related boundary conditions* The boundary conditions for this 
approximate superposition analysis are calculated using a simple strength of 
materials analysis. The relative simplicity of the resulting equations per- 
mits general conclusions to be made. For example, the equations show clearly 
the dependence of Gj and G It on delamination length. Hence, results for 
one delamination length could be used immediately to calculate Gj and Gtj 
for other delamination lengths. 

The superposition analysis requires several constants which are 
calculated with a linear finite element analysis. In reference 8 these con- 
stants were independent of delamination length, applied load, and initial 
imperfection. A normalization technique developed in this study make the 
required constants independent of the flexural and axial stiffness of the 
buckled region. The original superposition analysis accounted for mechanical 
loads, whereas the new version also includes thermal loads. 

The objectives of this paper are to describe the enhanced superposition 
analysis and to use this analysis to predict the effect of various parameters 
on instablity-related delamination growth. All conclusions are based on the 
assumption that G-j- and Gjj govern delamination growth. 
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NOMENCLATURE 


a 

Aa 

b 

D 


half-length of delamination before loading, m 

virtual crack closure distance used in strain-energy release rate, 
calculations, m 

specimen width, m 

bending stiffness of the buckled region given by 


D 


t, P 

I l 

k-1 


(*k ~ 2 ) " (Vl 



where p - number of plies 


d xl» d yl unit load solutions for displacements near crack tip for load case 
(P c - P D ), m/N 

d v ? ,dy? unit load solutions for displacements near crack tip for load case 
M, IT 1 


d xl ,d yl 

normalized 

values of 

d xl and 

d yl , mN“ 1/2 

d x2 ,d y2 

normalized 

values of 

d x2 and 

d y2 , mN“ 1/2 

E k 

Young's modulus for ply k, MPa 



F 11 ,F 22 ,F 33 Young's moduli of unidirectional ply. The subscripts 1, 2, 
and 3 refer to the longitudinal, transverse, and thickness 
directions, respectively, MPa 


F xl> F yl 
F x2’ F y2 
F xl ,F yl 


F x 2 ,F y2 


unit load solutions for forces at crack tip for load case (Pq - Pjj) 
unit load solutions for forces at crack tip for load case M, m - ^ 

„l/2 
1/2 


normalized values of F x ^ , , 


normalized valules of F x 2 , F y 2 , N 
mode I strain-energy release rate, j/m 2 
mode II strain-energy release rate, j/m 2 
®12>®13»®23 shear moduli for unidirectional ply, MPa 

Gjj strain-energy release rate related to moment at crack tip, J/m 2 


'II 


Vl 


total strain-energy release rate, J/m 2 

distance from top surface of laminate to ply "k"; top ply is 
ply 1, m 
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M 


moment , Nm 


P A> P B’ P C> P D axial loads carried in regions A, B, C, and D, respectively, N 

P T remote applied compressive load, N 

mechanical load for bifurcation buckling when there is no thermal 
load, N 

S A» S B» S C> S D axial stiffness of regions A, B, C, and D given by 


s " b I Ek (\ - Vi) 

k=l xxi 
where p = number of plies , N 

AT increase in temperature, °K 

t B ,t c thickness of regions B and C, respectively, m 

x,y rectangular Cartesian coordinates, m 

a A» a D axial coefficients of thermal contraction for regions A and D, 
respectively, °K -1 

8 thermal load, N 

6 lateral deflection at x - -a due to applied load, m 

6 q lateral deflection at x =■ -a before application of mechanical or 

thermal load, m 


axial strain; compressive strain is defined to be positive 


ANALYSIS 

Two analyses were used in this study: geometrically nonlinear, two- 

dimensional finite element analysis and an approximate superposition 
analysis. A geometrically nonlinear finite element analysis was used to 
illustrate the sensitivity of Gj and G IX to various geometric and material 
parameters. It was also used as a check for a simpler, but less rigorous, 
superposition analysis. Details of the geometrically nonlinear analysis are 
given in reference 7. The finite element models were similar to those 
described in reference 8. 
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The second analysis used was an approximate superposition analysis. This 
analysis is an extension of that presented in reference 8. In the following 
subsections, the approximate superposition analysis will be discussed first. 
Then the procedure for calculating Gj and Gjj is described. Then the 
normalized constants for the superposition analysis are derived. Normalized 
strain-energy release rate curves for Gj and G^-j- versus applied load are 
described for the case of no initial imperfection. Finally, the material 
properties are discussed. 

Approximate Superposition Analysis 

Superposition techniques have been widely used in linear stress analysis 
to represent a complicated problem as a combination of several simpler prob- 
lems. Application of the principle of superposition to nonlinear problems 
first requires a transformation that results in a linear system. 

The key to the transformation is replacement of the source of nonlinear- 
ity with equivalent loads (figs. 2(a) and (b)). Because of symmetry only half 
of the configuration is considered. The buckled region (which responds non- 
linearly due to significant rotations) is replaced by the loads Pg and M, 
the axial load and moment, respectively, in the column where it is cut 
(fig. 2(b)). The governing equations for the new conf iguration are linear. 
There are four nonlinear ly related loads: three mechanical loads P T , P D , 

and M, and internal loads due to the thermal initial strain. Figure 2(b) 
shows the mechanical load P T replaced by Pg and P c , the portion of the 
applied load carried by regions B and C, respectively. The magnitudes of 
Pg and P c are determined by solving equations (la) through (lc). 

P B + P C = P T 
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(lb) 


P B = S A (£ " «A AT) 

P c = S D (e - AT) (lc) 

Equations (1) are based on the assumption that the axial strain e is 
constant through the thickness of regions B and C. This assumption is valid 

if there is negligible bending in regions B and C. 

The load system in figure 2(c) (which is the same as in fig. 2(b)) can be 
divided into the two load systems shown in figures 2(d) and 2(e). Because 

Pq and P s are calculated using rule of mixtures (including thermal 

effects) , the load system in figure 2(e) causes a uniform axial strain state 
and no interlaminar stresses. Consequently, to calculate interlaminar 
stresses, only the load system in figure 2(d) (i.e., (Pq - Pp) and M) need be 
considered. Accordingly, in the current study involving strain-energy release 
rates, the configuration in figure 2(d) is the linearized equivalent of the 
nonlinear problem in figure 2(a). Note that the configuration in figure 2(d) 
has no initial strains due to thermal loads; all relevant thermal effects are 
included in Pq. 

Appendix A describes a strength of materials analysis for calculating 
Pq, P d , and M. The key equations from the appendix are 


tt 2 S a (A 2 + 266 0 ) D(S a + S D ) 


T 2 16 

a 


S A S D 5 + 5 o 


‘ ( “A ‘ °D> 


ATS. 


P C = S. + S t P T + S A AT( a A ~ 
A D 


tt 2 D 6 


P D 26+6 
a o 


M = 


tt 2 D 6 
2a 2 


( 2 ) 

(3) 

(4) 

(5) 
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Equations (2) through (4) can be combined to obtain an expression for P q - Pp 



P - P = 


16a 2 S A 


To use the loads (Pq - p d> and M in a two-dimensional analysis requires that 
they be expressed as an equivalent distribution of tractions. To calculate 
this distribution, the axial strains were assumed to vary linearly through the 
thickness where the tractions are applied (i.e., at the cut). Intuitively, 
this seems to be reasonable if region D (fig. 2) is not cut too close to the 
crack tip. The validity of the assumed linear variation was established in 
reference 8. 


Linear finite element analysis was used to calculate the response of the 
linearized configuration in figure 2(d) to unit values of (Pq - P^) and M. 
Because the configuration is linear, the solution for any arbitrary combina- 
tion of (Pq - P D ) and M is simply a linear combination of the unit load 
responses. If region B (fig. 2) is much thicker than region C, the unit load 
solutions are very insensitive to delamination length. In the current study 
the ratio of thicknesses, tg/tQ, was always large (i.e., at least 9/1). As 
demonstrated in a later section, the unit load solutions can be normalized in 
such a way that the normalized parameters become independent of the in-plane 
and flexural stiffness of the buckled region. 

Delamination length, initial imperfection, applied load, and temperature 
changes are all accounted for in the strength of materials analysis in calcul- 
ating (Pq - P D ) and M. 
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Strain-Energy Release Rate 

The virtual crack closure method (ref. 16) was used to calculate mode I 
and mode II strain-energy release rates, Gj and Gjj, respectively. The 
forces transmitted through the node at the crack tip and the relative dis- 
placements of the two nodes on the crack boundary closest to the crack tip 
were used in the calculation. Equations (7) show how this technique is used 
for the superposition stress analysis. 

°I - lW (P C - V F yl + «V][<V W + “Vl 

( 7 ) 

°II - ib [< P C - V F xl + ] [< P C - V d *l + Md x2 ] 

In these equations F x , Fy, d^, and dy are the unit load values of the 
nodal forces and the corresponding relative nodal displacements in the 
x— and y-directions. (The coordinate system is defined in fig. 2.) 

The superscripts 1 and 2 on the unit load parameters identify parameters 
associated with (P c - P D ) and M, respectively. 

If the distance is small between the crack tip and the nodes used to 
calculate relative displacements (which it should be for accurate results) , 
then F j/d ^ = F^/d^ and F xl^ d xl " F x2^ d x2* Usin 8 the se relationships in 
equations (7) results in 

G I = ISEf^[< P C- P D> F yl + MF y2] 2 

( 8 ) 

G XI ■ V F xl +MF x 2 ] 2 

xl 
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Reference 8 gives the unit load solutions for one laminate. The next section 

describes a normalization procedure which allows use of these same values for 

\ 

many different laminates. 


Normalized Unit Load Solutions 

The effect of the buckled region's in-plane and flexural stiffness on the 
unit load solutions can be approximated using strength of materials analysis. 
This approximation depends on the following two assumptions: 

1. The crack tip compliances (d^/F^, d j/P ^ etc.) do not vary with 

the in-plane and flexural stiffnesses of the buckled region. Hence, 
the material properties immediately around the crack tip should not 
be changed. 

2. The ratios G i/ G u for ( p c “ P]j) and M applied individually are 

independent of the in-plane and flexural stiffnesses of the buckled 
region. 

Implicit in these assumptions is that region A is much thicker than region D. 
The accuracy of these assumptions will be evaluated later. The required 
normalization is different for the two load cases (P^ - P D ) and M, so each 
load case will be examined individually. 


Load case 1: (P c - P Q ) 

Using simple strength of materials analysis (see ref. 17) one can show 
that the total strain-energy release rate for the load case (Pq - P D ) is 


G 


T 



2b 


S A + S D 
S A S D 


(9) 
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From equations (8), the total strain-energy release rate is 


■ °I + G II ■ 2k ( p c - p D ) 2 


zL fF f + f P f 

F , yl J F . ^ xl J 


yi 


xl 


( 10 ) 


Combining equations (9) and (10) yields 


2 Aab 


yi 


xl 


( F yl ) + FT ( F xl ) 


F yl ^ yl "xl 


S A + S D 
S A S D 


(ID 


From assumption (1), the compliances d ^/F^ and d^/F^ are constant. 
Since these compliances are constant, assumption (2) requires that 


F ^ and F x ^ change by the same scale factor as S A and S D are changed. 

Hence, the normalized unit load solutions F , and F , , which are 

yl xl 


independent of S A and S D , are 

/ Va f 2 


5,1 yl \ S A +S D^ 


and F 


S A S D 


.1/2 


xl xl \ S . + S_ 
A D y 


( 12 ) 


Since the crack-tip compliances d^/F^ and d^/F^ are assumed to be 
constant, the normalized unit load displacements are 


d . - d . 
xl xl 


S A S D 


1/2 


S A +S D 


and d , = d 


S A S D 


yl/2 


y! y 1 \ S A + S D 


(13) 


Therefore, if the unit load solutions are known for one combination of S A 
and Sq, the normalized unit load forces and displacements can be calculated 
using equations (12) and (13). Then these normalized parameters can be used 
with equations (14) to calculate the unit load solutions for other combina- 
tions of S A and Sjj. 
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1/2 


F - F 


S A + S' 


Xl Xl \ S A S D 


d , = d 


<S A + S D 


1/2 


Xl Xl \ S A S D 


F „ - F 


' S A + S 


1/2 


yi yi\ s A s D 


d_ , - d 


S A + S' 


1/2 


y! yi \ s A s D 


( 14 ) 


Load case 2: M 

Again, simple strength of materials arguments are used to derive an 
expression for the total strain-energy release rate (eq. (15)) 



(15) 


Following the same procedures used for load case 1, the normalized unit load 
solutions are found to be 






V 


1/2 




(16) 


After the unit load solutions are calculated for one value of D, the normal- 
ized unit load solutions can be calculated using equations (16). These 
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normalized solutions can be used with equations (17) to calculate the unit 
load solutions for other values of D. 


' x2 


F x2 D 


- 1/2 


x2 


y2 


d x2 D 


- 1/2 


(17) 


V ° 


- 1/2 


y2 


d y2 D 


- 1/2 


Therefore, the mode I and mode II strain-energy release rates for a range 
of S A , S D , and D can be calculated by combining equations (8), (14), and 
(17). 


1 d . 

Q SB — X 

I 2 Aab F . 

yi 


v <&r 


+ MF . D 
y2 


1/2 


(18) 


xl 


7 II 2 Aab F 


xl 


S * + s d\ 


1/2 


(P c" P d )f xi\ s a s d j + 


mf x2 D 


1/2 


The accuracy of these derived unit load solutions depends on how much 
different the current configuration is from the configuration which was ana- 
lyzed with the finite element analysis. Obviously, one should use the finite 
element analysis on the configuration which lies in the middle of the range of 
interest. The unit load solutions for the laminates considered herein are 
given in Table 1. 

Note that equations (9) and (15) can be summed to obtain a very simple 
expression for the total strain-energy release rate, G-j. for combined loading 



( 19 ) 



This simple summation is possible because there is no coupling between 
(P c - Pp) and M in terms of the potential energy of the applied loads. 

Normalized Strain-Energy Release Rate Curves for the Case 6 Q = 0 
In the preceding sections equations are presented which allow calculation 
of Gj and G II* These equations are easily solved for any applied load 
(mechanical or thermal or both). Additional insight is possible if parameters 
are identified which coalesce the G^ and G^j vs. applied load curves for 
various configurations into a single curve. Reference 8 showed that by plot- 
ting vs. a^Gj, the Gj vs. P^> curves for various delamination 

lengths coalesce into a single curve. In the current study parameters were 
derived which provide a single curve for arbitrary length and in-plane and 
flexural stiffnesses of the buckled region. Also, the loading can be 
mechanical, thermal, or combined mechanical and thermal. The curve derived 
herein is only valid when the initial imperfection 6 Q is zero. 

The key to deriving the normalized curves is expressing P T , (P c - P Q ) , 
and M in terms of the normalized load parameter P T (eq. (20)): 

P = _iL i (20) 

T P C 
TM 

where g - (a A - a^) ATS^ - thermal load 
P T => mechanical load 

P,^ = mechanical load for buckling for case in which g = 0 
The parameter P^ M is determined from equation (1) by substituting 
6 “ 6 q * AT ■ 0. The result is 
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( 21 ) 



2_ S a + S n 

it D A D 


Combining equations (2), (20), and (21) yields an expression relating 
normalized load to lateral deflection. 


■ - 6 2 

T 16D(S a + S D ) 


( 22 ) 


Equations (3), (4), (21), and (22) can be combined to yield an expression for 
(P C - P D ) in terms of P^. 


2 

(p _ p ) = jl 2. p 
' C 2 T 

a 


(23) 


Equations (5) and (22) can be combined to yield an expression for M in terras 

m/2 


of P_ . 
T 


M = 


2tt 2 D 


D(S A + V 


S A S D 


pi /2 

T 


(24) 


If the expressions for (P c - Pp) and M in equations (23) and (24) are com- 
bined with equation (18), we obtain equations for Gj and Gjj in terms 
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( 25 ) 


; -» — 1 d yl Tt 4 d 2 S A + S D /- - 2=1/2- V 

1 2Aab F yl a 4 S A S D ' T yl T y2 ' 


G ■ -i— ^l^D 2 !A_ + iD /g 
11 2 AabF xla 4 S A S D V 


F + 2P^ 2 F ] 
T xl T x2J 


D (S A + V 

Next, both sides of equations (25) are divided’ by r to obtain 

a 4 S.-S_ 

AD _ 

normalized strain-energy release rate parameters which are functions of P 
alone. 


G I a 


-JjLR iLik (p F + 2 P 1/2 F V 

D (S A + V 2Aat Fyl T ^ T ^ 


a 


4 S A S D 


iL^Wp 5 + 2P 1/2 F ) 2 

' IX “ D 2 (S + SJ 2Aab F xl \ T Xl T 


'A D 


(26) 


3 ^A^D 3 ^A^D - 

Hence, the G_ — = or G TT — = vs. P curves for different 


°< S A + V 


D <S A + S D ) 


configurations should all coincide. This will be verified in the Results and 
Discussion section. 


Materials Properties 

The unidirectional ply properties were assumed to be 
E n - 140 GPa 
E 22 3 E 33 * 14 GPa 
v 12 “ v 13 " v 23 " °* 21 
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G ^2 * G ^3 = G 23 — 5.9 GPfl 

Plane strain (i.e., e z = 0) and * 0 were imposed to calculate the 2D 

properties. In regions where coarse finite elements spanned several plies, 
laminate theory was used to obtain average properties. Lamina thickness was 
assumed to be 0.14 mm. 

RESULTS AND DISCUSSION 

A parametric study of instability-related delamination was performed. 

The parameters considered are delamination length, in-plane and flexural 
stiffness of the buckled region, initial imperfections, mechanical and thermal 
load, and interlaminar fracture toughness. Four laminate types are considered. 
These are described in Table 1. All conclusions are based on the assumption 
that Gj and Gjj govern delamination growth. 

First, Gj and Gj-j will be discussed for a variety of cases. When 
possible, simple expressions will be presented to quantify the effect of 
various parameters on Gj and G^j. When the effect of a parameter is not 
described by simple expressions, example results will be discussed which 
qualitatively illustrate the trends. 

Next, the effect of various parameters on delamination growth is examined 
analytically. Several delamination growth criteria are considered. The 
assumed values of G Ic and Gj Ic range from values representative of brittle 
epoxy systems to tough thermoplastic systems. 

Parametric Study of Gj and Gjj 

When the initial imperfection 6 Q is zero, the equations for P^, Pq, 
and Pjj (eqs. (2), (3), and (4)) simplify considerably. The result is that 
simple equations for Gj and G^-j- (eqs. (26)) can be derived. These equa- 
tions show clearly the effect of various parameters on the relationship 
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between Gj, Gj-j, an< * normalized load P^,. When 6 q * 0, the expressions 
for Gj- and Gjj are not so simple. Hence, in this section the cases 
6 Q = 0 and 6 Q * 0 will be discussed separately. The case 6 Q = 0 will be 
discussed first. 

When 6 q = 0, equations (26) express the relationships between Gj, Gjj, 

and normalized load for arbitrary values of a, S^, S D , and D. For a 

° 2 (S A + S D ) 

given value of P^, the magnitudes of Gj and G^-j vary as ^ * 

a S A S D 

Hence, small changes in flexural stiffness D or delamination length 2a 
should result in large changes in Gj and Gjj. Changes in and Sp 

also should affect Gj- and G n * Figure 3 verifies this prediction. 

Figure 3(a) shows the effect of delamination length on G^ and Gjj and 

figure 3(b) shows the effect of changing D, S^, and Sp through variations 

in laminate type. These results were obtained with geometrically nonlinear 
finite element analysis. The results are plotted on semi-log scales — hence, 
the differences between the curves represent very large differences. Figure 3 
shows that Gj first increases, then decreases with increasing applied 
load. This peculiar behavior is explained in reference 7. Figure 4 shows 
that equations (26) accurately quantify the effects of D, a, S^, and Sp 


S A S D 


on Gj and Gjj. By plotting a — = G T 

d ( s a + s d) 


vs . P^, and 


Vd 


D ( S A + S d) 


Gjj. vs. P^, all of the data in figures 3(a) and 3(b) have 


been reduced to one curve for Gj and one for G XI . The curves in figure 4 
are plotted on linear scales. Hence, the scatter is quite small. Also, 
figure 4 shows that if nonlinear finite element results are available for one 
configuration, the Gj- and G^j values for other configurations can be 
estimated immediately. 
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Equation (2) shows that the thermal load f} 31 (a A - a D )ATS A simply 
shifts the relationship between Fj and <S. This shift causes a shift in 
the Gj and Gjj vs. relationship. Figure 5(a) shows finite element 

results for two different thermal loads. The two curves have nearly the same 
shape, but they are offset by the magnitude of the thermal load. The normal- 
ized load P T includes the effect of the thermal load, hence the two curves 
agree quite well when normalized Gj vs. ? T is plotted (fig. 5(b)). In 
figure 5(a), for one of the cases G-j- is nonzero when P T » 0. This simply 
means that the thermal load alone was sufficient to buckle the sublaminate. 

In fact, for large thermal loads, delamination growth may occur without any 
applied mechanical load. 

As previously noted, when * 0 the expressions for Gj and Gjj are 
not nearly as simple as when 6 Q » 0. Consequently, the effect of 6 Q on 
Gj- and Gjj; is not obvious from the govenlng equations. Using figures 6 
through 10 the effect of S Q on G^ and G^-j is qualitatively assessed by 
examining particular results of the approximate superposition analysis. 

Herein, the term "initial imperfection" refers to the lateral distortion 
which would exist if the specimen was stress-free, i.e., without mechanical, 
thermal, or hygroscopic stresses. A possible cause of an initial imperfection 
is an inclusion left between lamina during fabrication. 

Figure 6 shows the effect of S Q on Gj and G-j-j. Because normalized 
strain-energy release rates are plotted vs. normalized load, the curves for 
different delamination lengths are the same. There are, however, different 
curves for different values of 6 Q . Note that Gj is very sensitive to 6 Q , 
whereas G-j-j is almost insensitive to 6 Q . The peak value of Gj decreases 
rapidly as 6 Q increases. In fact, if P^, is greater than about 0.2, Gj 
decreases with Increased 6 Q . When P,j, is negative, the applied load is 



insufficient to cause bifurcation buckling, i.e., buckling of a perfectly 
straight column. But when <S q * 0, lateral deflections occur as soon as 
compressive stresses develop in the delaminated region. Concomitantly, non- 
zero values of G-j- and occur for -1 < P^, < 0. In this region Gj 

and Gjj may increase with increased 5 Q . Therefore, if a delamination is 
short and is located deep in the interior of a laminate and if <S Q “ 0, then 
the operating strain is unlikely to be sufficiently large to cause bifurcation 
buckling. Hence, G^ and G II would be zero and delamination growth would 
not occur. But if * 0, Gj and Gj-j- can exist even for relatively small 
strain, and growth might occur. 

Figure 7 shows the effect of S Q on Gj and Gjj for three laminates 
with different values of S A , S D , and D. In figure 7 the solid curve is for 
the case 6 Q = 0. The dashed curves are for the case 6 Q = 0.1 mm. As shown 
earlier, the curves for all three configurations are the same when <$ 0 ■ 0. 
When S Q * 0 there is a different curve for each combination of S A , S D , 
and D. The effect on Gj is much larger than the effect on G^j. As the 
flexural stiffness of the buckled region Increases the effect of 6 Q on Gj 
decreases. (The flexural stiffness is largest for laminate type 4 and 
smallest for laminate type 1.) 

Figure 8 shows the effect of 6 Q on the Gj/Gjj ratio. For all 


values 

of 

^T 

the 

ratio 

Gi/Gii 

decreases as 6 Q increases. For all 

values 

of 

6 o 

the 

ratio 

Gj/Gh 

decreases rapidly as P T increases. There 


fore, whether delamination growth is governed by Gj or Gjj is likely 
dependent on the magnitudes of 6 Q and P T . 

Earlier it was pointed out that thermal loads can cause buckling when the 
applied load is zero. The lateral deflection due to thermal load may be 

visually indistinguishable from an initial imperfection 6 Q . But the thermal 
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load $ affects Gj and Gjj differently than 6 Q . Figure 9 illustrates the 
difference. For reference the case 6 0 = 0, 6=0 is also included. The 

magnitudes of 6 Q and 3 were chosen such that when P-j, = 0, the lateral 

distortion was the same for both cases. The thermal load case has a larger 

value of Gj than the initial imperfection case for all values of P^,. The 

value of Gjj was also found to be larger for the thermal case. (This result 

is not included in fig. 9.) The lateral deflection increases more quickly for 
the thermal case than for the intial imperfection case. 


Prediction of Delamination Growth 

To predict the critical load for delamination growth (herein referred to 
as the "critical load”) requires three types of information: (1) the 

relationship between applied load and strain-energy release rates, (2) the 


growth criterion (e.g., 1^ 


Ic, 


II 

; ilc- 


=11, and (3) the material constants 


to be used in the growth criterion (e.g., Gj c and Gjj c ). The relationship 
between applied load and strain-energy release rates can be determined using 
the analysis described herein. At this time there is no widely accepted 
delamination growth criterion for mixed- mode configurations. Hence, several 
criteria will be considered. The three criteria selected are given in equa- 
tions (27) through (29). 
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( 27 ) 


(28) 
(29) 

The criterion = Gj c assumes that only the mode I component is impor- 
tant. The criterion G^ » Gj c also assumes that only the mode I component is 
important. The quantity G^ is the strain-energy release rate due to the 
moment M at the delaraination front. This criterion was proposed in refer- 
ence 12. The mixed-mode criterion, equation (29), attempts to account for 
both the mode I and mode II effects. The mixed-mode criterion will be 
emphasized in the following discussion. 

The material constants in the growth criteria considered consist of 
only G Ic and Gjj c . Numerous studies have reported values of G Ic (e.g. , 
refs. 10, 15, 18, 19). These Gj c values range from about 200 J/m^ 
to 1000 J/m^. The mode II fracture toughness Gjj c is less well 
characterized, but is generally assumed to be several times as large as 
G Ic . Herein, the ratio G IIc /G Ic ran 8 e s from 1 to 5. The critical load is 
determined from the growth criterion equation and the equations for strain- 
energy release rates (i.e., Gj, G^j, or Gj^). The procedures for solving 
the equations are given in Appendix B. 

Figure 10 shows the variation of critical load with delamination length 
predicted by the three criteria in equations (27) through (29). In figure 11 
the initial imperfection <$ 0 = 0. The assumed values of Gj c and Gjj c were 
200 j/m^ and 1000 J/m^, respectively. The laminate is type 2, which has a 
[0 4 ] buckled sublaminate. For short delaminations , all three criteria 
initially predict a reduction and then an increase in critical load with 
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increase in delamination length 2a. But beyond about 22 mm, even the trends 
are very different for the three criteria. The criterion G I “ G lc predicts 
no growth, regardless of the applied load, beyond about 23 mm. The 
criterion Gj^ = Gj c predicts an increasing but finite critical load. The 
mixed-mode criterion predicts an almost constant critical load after the 
initial growth. The mixed-mode criterion and the criterion G^ = G Ic both 
predict that if a short delamination grows, there will be substantial unstable 
growth before arrest. The differences in the predictions from the three 
criteria highlight the need for a verified growth criterion. 

The bifurcation buckling curve is shown as a dashed line in figure 10. 

For short delaminations there is little difference between the buckling load 
and the critical load for delamination growth. From a practical standpoint, a 
short delamination must not be permitted to buckle — if it buckles, one must 
assume it will grow. 

Figure 11 shows results analogous to those in figure 10 except that 
6 Q = 0.1 mm. For short delaminations (i.e., < 15-20 mm), the initial imper- 
fection reduces the critical load from that shown in figure 10. For long 
delaminations the Initial imperfection has little effect. The criterion 
Gj = Gj c and the mixed-mode criterion agree for short delaminations (see 
figs. 10 and 11). Hence, the Gj component must dominate in the mixed mode 
criterion for short delaminations. For longer delaminations the Gjj compo- 
nent dominates. The bifurcation buckling curve from figure 10 is shown in 
figure 11 for reference. (Actually, bifurcation buckling does not occur when 
there is an initial imperfection.) For short delaminations the bifurcation 
buckling load is larger than the critical load for growth. Therefore, if 
6 Q ^ 0 it is unconservative to use the buckling load as an upper bound on the 
allowable load. 
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To expedite the remaining discussion, only the mixed-mode criterion will 
be used in the remainder of this paper. 

The curves in figures 10 and 11 depend on the particular specified values 
of Gj c and Gjj c . Of interest is how much the critical load increases for 
increases in G Ic and G IIc . Figure 12 shows the critical load predicted by 
the mixed-mode growth criterion for three combinations of G-j- c and G IIc . 

The initial imperfection 6 Q is zero for the cases in figure 12. For short 
de laminations large improvements in G Ic only moderately increase the critical 
load and Gjj c has no effect. For long delaiqinations, G^ c has no effect 
but increases in Gjj c result in substantial increase in the critical load. 
Figure 14 shows predicted critical loads for the case 6 Q = 0.1 mm. Compari- 
son of figures 12 and 13 shows that the initial imperfection only affects the 
predictions for short delaminations; for short delaminations the initial 
imperfection reduces the critical load. 

Uncertainty about the depth of a delamination (i.e., the thickness of the 
buckled region) and the magnitude of the initial imperfection causes signifi- 
cant uncertainty about the load at which a delamination will grow. Suppose 
that the buckled region has either 3 or 4 zero degree plies (laminates types 1 
and 2) . Figures 14 and 15 show the predicted critical load (using the mixed- 
mode criterion) for each possibility. Figure 14 shows predictions for the 
case ■ 0. For short delaminations the thicker buckled region has a higher 

critical load, but the trend is reversed for long de laminations. Figure 15 
shows results analogous to those in Figure 14 except that 6 Q ■ 0.1 mm. 
Comparison of figures 14 and 15 shows that the initial imperfection affects 
only the short delamination predictions. In the short delamination region 
the thicker laminate is affected more than the thinner one. In fact, for 
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6 Q = 0.1 mm the predicted critical loads are nearly the same for the two 
laminates in the short delamination region. 

The results presented here should be helpful in designing specimens to 
test various delamination growth criteria. The sensitivity of Gj and Gjj 
to various parameters illustrated in the figures indicates the accuracy with 
which these parameters must be known to obtain accurate values of G-j- and 

G II* 

CONCLUSIONS 

A parametric analytical study of instability-related delamination growth 
was conducted. All conclusions were based on the assumption that the mode I 
and mode II strain-energy release rates (G^ and Gjj* respectively) govern 
delamination growth. The configuration studied consisted of a thick composite 
laminate with a single delamination located near one surface. 

The primary accomplishments and conclusions from this study are: 

1 . An approximate superposition analysis was developed which gives 
closed form expressions for Gj and Gjj for a wide range of delamination 
length and location through the thickness, initial imperfection, and 
mechanical and thermal load. 

2. The magnitude of Gj is very sensitive to delamination length, 
initial imperfection, and delamination depth. For zero initial imperfection, 
simple expressions were derived which quantify the sensitivity to delamination 
length and delamination depth. 

3. The ratio Gj/Gjj decreases monotonically with applied load and 
initial imperfection. 

4. Initial deflection due to residual thermal stresses is worse than 
initial deflection due to initial imperfection. 
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5. Three growth criteria were considered. The large differences in the 
predicted critical loads highlights the need for a verified mixed-mode growth 
criterion. 

6. Short delaminations can be expected to grow if buckling occurs. This 
growth is likely to be rapid and extensive. 

7. Initial imperfection can cause a design based on "no bifurcation 
buckling” to be very unconservative. 

8. Based on the mixed-mode criterion: 

(a) A large increase in G^ c results in a moderate increase in 
critical load for delamination growth for short delaminations. 

(b) A large increase in results in a substantial increase in 

critical load for delamination growth for long de laminations. 
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APPENDIX A 


STRENGTH OF MATERIALS ANALYSIS 

A strength of materials analysis is described herein for the configura- 
tion in figure 1. 

The configuration is divided into four regions, as shown in figure 2. 
Because of symmetry, only half of the laminate is modeled. The laminate has 
width b. The following assumptions are made: 

1. Regions B and C are perfectly bonded. Regions A and D are unbonded. 

2. Regions A, B, and C have constant axial strain. Hence, the force- 

strain relations are those for a simple rod subjected to axial load 
and a temperatute change; e.g., 

i 

P A “ S A (e ' “A AT) (A1) 

Since regions A, B, and C are not permitted to bend, global bending 
effects cannot be predicted by this Idealization. 

3. Region D has zero slope at both ends. 

4. Region D is symmetric with respect to its midplane. 

5. Regions A and B have a thermal contraction coefficient a^. 

6. Regions C and D have a thermal contraction coefficient ajj. 

7. All regions experience the same AT. 

8. Region D has an initial sinusoidal imperfection of peak magnitude 

6 q . The term "initial imperfection" refers to the lateral 
distortion which would exist if the specimen was stress-free, i.e., 
without mechanical, thermal, or hygroscopic stresses. The initial 
shape is given by 
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v<x) I initial 'I^i 1 ' c ° s r) 


(A2) 


where v(x) =* the distortion in the y-direction. 

To describe the nonlinear behavior of region D, equations (A3) and (A4) 
for post-buckling of a column were used. 


n^D 6 

2 6+6 
a o 


(ref. 20) 


a - » - fer ( {2 + 255 <,) + + aa D 4T 


(A3) 


(A4) 


where 6, a, a, and P D are peak lateral deflection, axial length before 
and after deformation, and load, respectively. Compression strain is defined 
to be positive. The derivation of equation (A4) for mechanical loads only can 
be formed in reference 12. The additional tern for thermal loads, aopAT, 
requires no derivation. Equations (A3) and (A4) were derived using strength 
of materials analysis of a column. 

To combine regions A, B, C, and D, equilibrium and compatibility condi- 
tions must be considered. The equilibrium condition for the axial force is 


P A + P D = P B + P C " P T 


(A5) 


Compatibility requires the shortening of regions A and D to be identical. 
Hence, 

P A a 

a. ATa + •= — - a - a 
A S A 


(A6) 


Equations (A3) through (A6) can be combined to obtain the governing equation 
for the laminate. 
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2 C 

if s. 


P„ = 


; , d(s. + s_) 

i F ( 6 2 + 266 j+ A " 


S A S D 


6+6 


- (« A - c^) ATS a (A7) 


For a specified load P T , equation (A6) is solved using a Newton-Raphson 
technique to obtain 6. Pp can then be calculated using equation (A3). From 
static equilibrium, the moment acting on the delaminated region at the crack 
tip is 

P„ 


2 

it D 


M ._ (5+v ._ 


(A8) 


The force P^ is found from rule of mixtures by solving equations (A9): 


P B + P c = P T 


P B " S A< e ~ «A AT > 


(A9) 


P C “ S D^ e ~ a D AT ^ 


where e is the axial strain in regions B and C. The result is 


p c - STTT- [r T - h - «A>] 

A D 


(A10) 
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APPENDIX B 


CALCULATION OF CRITICAL LOADS 

This appendix describes the procedures used to calculate critical loads 
for delamination growth based on the three criteria in equations (27) through 
(29). The procedure is different for each criterion, so they are discussed 
separately. 


Criterion: G^ - Gj c 

By combining equations (5), (6), and (18), the mode I strain-energy 
release rate can be expressed as a function of 6. 

G I “ C l{ A l &2 + A 2 6 ) 


(Bl) 


where 


4 d . 

it yl 

4 F 

8a bAa yl 


S A S D 


8 (S a + S D ) yl 


A 2 - 2A X 6 q + DF y2 


For Gj =* Gj c equation (Bl) can be explicitly solved for 6. 



(B2) 


The critical lateral deflection is the smallest of the two values from equa- 
tion (B2). The corresponding critical load is determined by using the criti- 
cal 6 in equation (2). 
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Criterion: G M = G Ic 

By combining equations (5) and (15), the strain-energy release rate due 
to the moment can be expressed in terms of 6. 


" 8a 4 b 


(B3) 


For Gj^ * Gj c equation (B3) can be explicitly solved for 6. 


, _ 2a 2b 

6 "2“ x/d - G Ic 

it ’ 


The corresponding critical load is determined by using this critical S in 
equation (2). 

Criterion: g-) ♦ (^-) - 1 

The first step is to express both Gj and Gjj in terms of 6. The 
equation for Gj has already been derived - equation (Bl). The expression 
for Gjj is obtained in a similar manner and is given in equation (B4). 

.2 


G II “ C 2 ( *3^ + A 4 6 ) 


(B4) 


where 


Co = 


» U - 

IT Xl 

4 F 
8a bAa xl 


S A S D - 

8(S A + V Xl 


A 4 = 2A 3 6 q + DF x2 


Substitution of the expression for Gj (eq. (Bl)) and Gj-j- (eq. (B2)) 
into the failure criterion yields a nonlinear equation in terms of 6. 
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This nonlinear equation was solved using the Newton-Raphson method. An 
initial guess of 6 = 0.3 mm was used for much of the calculations. However, 
the stability of the iterative solution did not appear to be very sensitive to 
the initial guess for 6. Once 6 was determined, the critical load was 
calculated using equation (2). 
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TABLE 1.- LAMINATE PROPERTIES AND UNIT LOAD 
SOLUTIONS FOR SEVERAL LAMINATE TYPES 

(b = 25.4 mm, Aa = .0254 mm) 



I 

Laminate type^ j 

1* 

2 

3 

4 

S A , N 

14.59E6 

13.45E6 

13.45E6 

13.78E6 

S D , N 

1.475E6 

1.967E6 

1.087E6 

2.950E6 

D, N «m 2 

■ 

.02159 

.05118 

.04546 

.1727 

li xL 

.0936 

.0827 

.108 

.0689 

F x 2> m_1 

531. 

345. 

366. 

188. 

F y! 

.0261 

.0231 

.0301 

.0192 

F y2 , nT 1 

-252. 

-164. 

-174. 

-89.1 

d xl , m*N -1 

1.40 x 10 -10 

1.24 x 10 -10 

1.62 x 10 -10 

1.03 x KT 10 

d yl , m»N -1 

2.97 x 10 -10 

2.62 x 10" 10 

3.43 x 10 -10 

2.18 x 10 -10 


*Unit load solutions for laminate type 1 are from reference 8. 
Solutions for other laminates were obtained using the technique 
described in the text. 


Laminate ■ type 

Laminate 

1 

[0 3 /0/(0/45/90/-45) 7s /0 4 ] 

2 

[0 4 /(0/45/90/-45) 7s /0/45 2 /0] 

3 

[0/90 2 /0/(0/45/90/-45) 7s /0/45 2 /0] 

4 

[0 6 /0 2( 0/45/90/ - 45 )6s /0 3 /45/90 2 /45/ °1 

i 

The delamination location is indicated by the arrow. 

. . __ 
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Figure 1. - Local buckling of laminate with through-width delamination 

(6 0 = initial imperfection) 


Y 


LO 

CD 



WITH INITIAL STRAIN 
DUE TO AT 


WITH INITIAL STRAIN 
DUE TO AT 


(a) 


(b) 



WITH INITIAL STRAIN 
DUE TO AT 


(0 


I — 


a 



NO INITIAL STRAIN 



NO INTERLAMINAR STRESSES 


WITH INITIAL STRAIN 
DUE TO AT 


(e) 


Figure 2. - Nonlinear configuration (a) transformed into linear configuration 
(d) with two nonli nearly related loads, (Pc - Pd) and M. 










P T 

(a) Normalized Gj vs. normalized load 


Figure 4. - Normalized plots of Gj and Gjj vs. applied load. These 

plots include all of the data in Figure 3. (Curves are visual 
fit of data.) 




(a) Gj vs. applied load. 

Figure 5. - Effect of thermal load on Gj vs. applied load. The thermal 

load is due to an initial strain (a - a )aT = .002. 

A D 

(delamination length = 25 mm, initial imperfection, & 0 = 0; 
laminate type 2) 


41 


fVo 

D 2 (S A +S D ) 



3000 
2500 * 
2000 • 
1500 ■ 
1000 • 
500 • 



(b) Normalized Gj vs. normalized load 


Figure 5. - Concluded. 



(a) Effect on normalized Gj 


Figure 6. - Effect of initial imperfection on normalized Gj and Gjj vs. 

applied load, (laminate type = 1; arbitrary delamination length) 



Figure 6. - Concluded. 
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Figure 7. 
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(a) Effect on normalized Gj 


Effect of laminate type and initial imperfection on 
normalized Gj and G[j. (arbitrary delamination 


e 










DELAMINATION LENGTH, mm 


Figure 10. - Critical load for delamination growth using three growth 
criteria. Buckling curve shown for reference. (6j c = 20( 
J/m^; G i j c = 1000 J/m^; laminate type = 2; initial 
imperfection 6 n =0). 




DELAMINATION LENGTH, mm 


Figure 11. - Critical load for delamination growth using three growth 
criteria. Buckling curve shown for reference. (Gi c = 
200 J/m^; G j j c = 1000 J/m^; laminate type = 2; initial 
imperfection 6 0 = .1mm). 





























DELAMINATION LENGTH, mm 


Figure 15. - Effect of laminate type on critical load predicted by mixed 
mode criterion. (Gj c = 200 0/m^; Gj j c = 1000 J/m^; 
initial imperfection. 6 0 = .1mm). 
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